**************************************************************
*HURRICANES AND GAS GOUGING - TWO-WAY FIXED EFFECTS CONSIDERATIONS
**************************************************************
frame copy default reg_twfe, replace
frame change reg_twfe

**
*Creating indicator for 'never treated' stations
frame copy reg_twfe twfe_control, replace
frame change twfe_control
gcollapse (sum) pre_hur hur post_hur, by(station_id) 
gen never_treated=0
	replace never_treated=1 if pre_hur==0 & hur==0 & post_hur==0
keep station_id never_treated
save $data_clean\never_treated, replace
frame change reg_twfe
frame drop twfe_control

*Merging never treated indicator
merge m:1 station_id using $data_clean\never_treated

*Collecting event estimates
gen event=""
gen event_hur_pre=.
gen event_hur_pre_ub=.
gen event_hur_pre_lb=.
gen event_hur=.
gen event_hur_ub=.
gen event_hur_lb=.
gen event_hur_post=.
gen event_hur_post_ub=.
gen event_hur_post_lb=.
 
********************************
*BONNIE-CHARLIE 
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*               in the event window
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 
	replace hur_landfall_d0=1 if date==`=td(12aug2004)'	 //Bonnie Charlie Landfall (8/12/2004)
gen hur_window=.
replace hur_window=1 if hur_landfall_d0==1
gen event_day=.
replace event_day=0 if hur_window==1
*
*45-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	*Year check
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	replace event_day=`n' if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*45-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	replace event_day=`t' if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if BONCHAR_landfall==1 & date==`=td(12aug2004)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if BONCHAR_landfall==1 //Flag all stations impacted by Bonnie-Charlie
tab storm_name if treat==1 & hur_window==1 //Confirm no other storms affecting treatment stations	
	replace treat=0 if FRANCES_landfall==1 //Drop treatment stations that were also subsequently hit by Frances
	replace treat=0 if IVAN_landfall==1 //Drop treatment stations that were also subsequently hit by Frances
	replace treat=0 if JEANNE_landfall==1 //Drop treatment stations that were also subsequently hit by Frances
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
tab storm_name if est_sam==1 //confirm no other storms affecting treatment stations	
*
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="BONNIE-CHARLIE" in 1

lincom pre_hur_twfe
replace event_hur_pre=r(estimate) in 1
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 1
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 1

lincom hur_twfe
replace event_hur=r(estimate) in 1
replace event_hur_ub=r(estimate) + 1.96*r(se) in 1
replace event_hur_lb=r(estimate) - 1.96*r(se) in 1

lincom post_hur_twfe
replace event_hur_post=r(estimate) in 1
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 1
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 1

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk
save $data_clean/twfe_stack, replace

frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe event_day
*
*
********************************
*FRANCES
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*               in the event window
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 
	replace hur_landfall_d0=1 if date==`=td(06sep2004)'	 //Frances Landfall (9/6/2004)
gen hur_window=.
replace hur_window=1 if hur_landfall_d0==1
gen event_day=.
replace event_day=0 if hur_window==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	replace event_day=`n' if f`t'.hur_landfall_d0==1 & y==year	
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	replace event_day=`t' if l`t'.hur_landfall_d0==1 & y==year	
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if FRANCES_landfall==1 & date==`=td(06sep2004)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if FRANCES_landfall==1 //Flag all stations impacted by Frances
	tab storm_name if treat==1 & hur_window==1
	replace treat=0 if BONCHAR_landfall==1 //Drop treatment stations also hit previously by Bonnie/Charlie
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="FRANCES" in 2

lincom pre_hur
replace event_hur_pre=r(estimate) in 2
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 2
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 2

lincom hur
replace event_hur=r(estimate) in 2
replace event_hur_ub=r(estimate) + 1.96*r(se) in 2
replace event_hur_lb=r(estimate) - 1.96*r(se) in 2

lincom post_hur
replace event_hur_post=r(estimate) in 2
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 2
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 2

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
********************************
*IVAN
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*               in the event window
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 
	replace hur_landfall_d0=1 if date==`=td(16sep2004)'	 //Ivan Landfall (9/16/2004)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if IVAN_landfall==1 & date==`=td(16sep2004)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if IVAN_landfall==1
	tab storm_name if treat==1 & hur_window==1
	replace treat=0 if BONCHAR_landfall==1	
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="IVAN" in 3
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 3
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 3
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 3
*
lincom hur
replace event_hur=r(estimate) in 3
replace event_hur_ub=r(estimate) + 1.96*r(se) in 3
replace event_hur_lb=r(estimate) - 1.96*r(se) in 3
*
lincom post_hur
replace event_hur_post=r(estimate) in 3
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 3
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 3

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*JEANNE
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*               in the event window
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(26sep2004)'	 //Ivan Landfall (9/26/2004)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if JEANNE_landfall==1 & date==`=td(26sep2004)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if JEANNE_landfall==1
	tab storm_name if treat==1 & hur_window==1
	replace treat=0 if BONCHAR_landfall==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="JEANNE" in 4
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 4
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 4
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 4
*
lincom hur
replace event_hur=r(estimate) in 4
replace event_hur_ub=r(estimate) + 1.96*r(se) in 4
replace event_hur_lb=r(estimate) - 1.96*r(se) in 4
*
lincom post_hur
replace event_hur_post=r(estimate) in 4
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 4
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 4

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*ARLENE
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*               in the event window
* 	NOTES - STATIONS HIT BY ARLENE WERE ALL HIT BY DENNIS TOO. LIMIT POST SAMPLE 
*            TO 28 DAYS POST SO DON'T INCLUDE DIRECT DENNIS IMPACTS.
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(11jun2005)'	 //Arlene Landfall (6/11/2005)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if ARLENE_landfall==1 & date==`=td(11jun2005)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if ARLENE_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="ARLENE" in 5
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 5
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 5
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 5
*
lincom hur
replace event_hur=r(estimate) in 5
replace event_hur_ub=r(estimate) + 1.96*r(se) in 5
replace event_hur_lb=r(estimate) - 1.96*r(se) in 5
*
lincom post_hur
replace event_hur_post=r(estimate) in 5
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 5
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 5

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*DENNIS
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(10jul2005)'	 //Dennis Landfall (7/10/2005)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if DENNIS_landfall==1 & date==`=td(10jul2005)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if DENNIS_landfall==1
	tab storm_name if treat==1 & hur_window==1
	replace treat=0 if ARLENE_landfall==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="DENNIS" in 6
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 6
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 6
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 6
*
lincom hur
replace event_hur=r(estimate) in 6
replace event_hur_ub=r(estimate) + 1.96*r(se) in 6
replace event_hur_lb=r(estimate) - 1.96*r(se) in 6
*
lincom post_hur
replace event_hur_post=r(estimate) in 6
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 6
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 6

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*KATRINA (FL)
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(25aug2005)'	 //Katrina (FL) Landfall (8/25/2005)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if KATRINA_FL_landfall==1 & date==`=td(25aug2005)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if KATRINA_FL_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="KATRINA (FL)" in 7
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 7
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 7
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 7
*
lincom hur
replace event_hur=r(estimate) in 7
replace event_hur_ub=r(estimate) + 1.96*r(se) in 7
replace event_hur_lb=r(estimate) - 1.96*r(se) in 7
*
lincom post_hur
replace event_hur_post=r(estimate) in 7
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 7
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 7

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*KATRINA (LA)
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(29aug2005)'	 //Katrina (LA) Landfall (8/25/2005)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if KATRINA_LA_landfall==1 & date==`=td(29aug2005)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if KATRINA_LA_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="KATRINA (LA)" in 8
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 8
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 8
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 8
*
lincom hur
replace event_hur=r(estimate) in 8
replace event_hur_ub=r(estimate) + 1.96*r(se) in 8
replace event_hur_lb=r(estimate) - 1.96*r(se) in 8
*
lincom post_hur
replace event_hur_post=r(estimate) in 8
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 8
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 8

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*RITA
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(24sep2005)'	 //Rita Landfall (9/24/2005)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if RITA_landfall==1 & date==`=td(24sep2005)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if RITA_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="RITA" in 9
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 9
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 9
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 9
*
lincom hur
replace event_hur=r(estimate) in 9
replace event_hur_ub=r(estimate) + 1.96*r(se) in 9
replace event_hur_lb=r(estimate) - 1.96*r(se) in 9
*
lincom post_hur
replace event_hur_post=r(estimate) in 9
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 9
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 9

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*WILMA
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(24oct2005)'	 //Wilma Landfall (10/24/2005)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if WILMA_landfall==1 & date==`=td(24oct2005)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if WILMA_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="WILMA" in 10
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 10
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 10
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 10
*
lincom hur
replace event_hur=r(estimate) in 10
replace event_hur_ub=r(estimate) + 1.96*r(se) in 10
replace event_hur_lb=r(estimate) - 1.96*r(se) in 10
*
lincom post_hur
replace event_hur_post=r(estimate) in 10
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 10
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 10

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*ALBERTO
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(13jun2006)'	 //Alberto Landfall (6/13/2006)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if ALBERTO_landfall==1 & date==`=td(13jun2006)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if ALBERTO_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="ALBERTO" in 11
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 11
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 11
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 11
*
lincom hur
replace event_hur=r(estimate) in 11
replace event_hur_ub=r(estimate) + 1.96*r(se) in 11
replace event_hur_lb=r(estimate) - 1.96*r(se) in 11
*
lincom post_hur
replace event_hur_post=r(estimate) in 11
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 11
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 11

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*HUMBERTO
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(13sep2007)' //Humberto Landfall (9/13/2007)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if HUMBERTO_landfall==1 & date==`=td(13sep2007)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if HUMBERTO_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="HUMBERTO" in 12
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 12
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 12
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 12
*
lincom hur
replace event_hur=r(estimate) in 12
replace event_hur_ub=r(estimate) + 1.96*r(se) in 12
replace event_hur_lb=r(estimate) - 1.96*r(se) in 12
*
lincom post_hur
replace event_hur_post=r(estimate) in 12
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 12
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 12

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*GUSTAV
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(01sep2008)' //Gustav Landfall (9/1/2008)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if GUSTAV_landfall==1 & date==`=td(01sep2008)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if GUSTAV_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="GUSTAV" in 13
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 13
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 13
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 13
*
lincom hur
replace event_hur=r(estimate) in 13
replace event_hur_ub=r(estimate) + 1.96*r(se) in 13
replace event_hur_lb=r(estimate) - 1.96*r(se) in 13
*
lincom post_hur
replace event_hur_post=r(estimate) in 13
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 13
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 13

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
*
********************************
*IKE
* Steps: (i) Define event window (30 days prior, 14 days after)
*        (ii) Define control (never treated stations)
*        (iii) Narrow treatment to stations not affected by another storm
*
*Hurricane landfalls day 0
gsort station_id date
gen hur_landfall_d0=0 	
	replace hur_landfall_d0=1 if date==`=td(13sep2008)' //Gustav Landfall (9/1/2008)
gen hur_window=0
replace hur_window=1 if hur_landfall_d0==1
*
*30-day Pre 
forvalues t = 1/45 {
	local n=-1*(`t')
	gen y=f`t'.year
	replace hur_window = 1 if f`t'.hur_landfall_d0==1 & y==year
	drop y
}
*
*30-day Post 
forvalues t = 1/45 {
	*Year check
	gen y=l`t'.year
	replace hur_window = 1 if l`t'.hur_landfall_d0==1 & y==year
	drop y		
}
*	
*Creating Hurricane, Pre-Hurricane, and Post-Hurricane Indiators
gen hur_twfe=0 
	replace hur_twfe=1 if IKE_landfall==1 & date==`=td(13sep2008)'
gen post_hur_twfe=0
gen pre_hur_twfe=0

*Filling in pre-hurricane indicators
gsort station_id date
gen hur_temp=0 //First day of hurricane
	replace hur_temp=1 if hur_twfe==1 & l.hur_twfe==0
forvalues t = 1/3 {
	gen y=f`t'.year
	replace hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=f`t'.year
	replace pre_hur_twfe=1 if f`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp
*
*Filling in post-hurricane indicators
gen hur_temp=0 //Last day of hurricane
	replace hur_temp=1 if hur_twfe==1 & f.hur_twfe==0	
forvalues t = 1/3 {
	gen y=l`t'.year
	replace hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
forvalues t = 4/14 {
	gen y=l`t'.year
	replace post_hur_twfe=1 if l`t'.hur_temp==1 & y==year
	drop y
}
*
drop hur_temp	
*
*Treatment and estimation sample
gen treat=0
	replace treat=1 if IKE_landfall==1
	tab storm_name if treat==1 & hur_window==1
gen est_sam=0
	replace est_sam=1 if (treat==1|never_treated==1) & hur_window==1
	tab storm_name if est_sam==1 
*
*Regression	
reghdfe retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN temp temp2 if est_sam==1, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
tab storm_name if e(sample)==1 //confirm no other storms affecting treatment stations

*Storing relevant estimates
replace event="IKE" in 14
*
lincom pre_hur
replace event_hur_pre=r(estimate) in 14
replace event_hur_pre_ub=r(estimate) + 1.96*r(se) in 14
replace event_hur_pre_lb=r(estimate) - 1.96*r(se) in 14
*
lincom hur
replace event_hur=r(estimate) in 14
replace event_hur_ub=r(estimate) + 1.96*r(se) in 14
replace event_hur_lb=r(estimate) - 1.96*r(se) in 14
*
lincom post_hur
replace event_hur_post=r(estimate) in 14
replace event_hur_post_ub=r(estimate) + 1.96*r(se) in 14
replace event_hur_post_lb=r(estimate) - 1.96*r(se) in 14

*Saving data for stacked dif-in-dif
frame copy reg_twfe temp
frame change temp
keep if e(sample)==1 
keep retail pre_hur_twfe hur_twfe post_hur_twfe wholesale CTST CHRT CTSN CHRN ///
	 temp temp2 est_sam station_id year month dow county_FIPS treat /// 
	 ALBERTO_landfall-WILMA_landfall never_treated state_name storm_name ///
	 nearRack1 bulk

append using $data_clean/twfe_stack
save $data_clean/twfe_stack, replace
frame change reg_twfe  
frame drop temp
drop treat hur_window est_sam hur_landfall_d0 hur_twfe post_hur_twfe pre_hur_twfe
*
*
********************************************************************************
*TABLE B.6 - Average Effect of Hurricanes on Retail and Wholesale Prices and Margins – Stacked Event-by-Event Estimates
********************************************************************************
frame create temp
frame change temp
use $data_clean/twfe_stack, clear

*****************
* PANEL 1a AVERAGE HURRICANE IMPACTS ON RETAIL AND WHOLESALE PRICES
*****************
eststo clear 		   
*Column 1 - Retail 1
reghdfe retail pre_hur hur post_hur CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(state_name year month dow) cluster(county_FIPS)	   
	gdistinct station_id if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "Yes"
	estadd local station_fe "No"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"
	estadd local dow_fe "Yes"
est store A
			   
*Column 2 - Retail 2
reghdfe retail pre_hur hur post_hur CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(station_id year month dow) cluster(county_FIPS)	
	gdistinct station_id if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "No"
	estadd local station_fe "Yes"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"	
	estadd local dow_fe "Yes"	
est store B

*Column 3 - Wholesale 1
	reghdfe wholesale pre_hur hur post_hur CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(state_name year month dow) cluster(nearRack1)
	gdistinct nearRack1 if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "Yes"
	estadd local station_fe "No"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"
	estadd local dow_fe "Yes"	
est store C
			   
*Column 4 - Wholesale 2
reghdfe wholesale pre_hur hur post_hur CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(nearRack1 year month dow) cluster(nearRack1)		
	gdistinct nearRack1 if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "No"
	estadd local station_fe "Yes"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"	
	estadd local dow_fe "Yes"	
est store D

esttab A B C D using "$output/twfe_a_$outputdate.csv", replace label ///
     b(a2) nonumber ///
    starlevels(* 0.10 ** 0.05 *** 0.01) ///
	title(Average Effect of Hurricanes on Retail Prices, Wholesale Prices, and Margins)  ///
	cells(b(fmt(3) star) se(fmt(3) par)) ///
	drop(_cons CTST CHRT CTSN CHRN temp temp2) ///
	note(Notes: The dependent variable is station-level retail/wholesale ///
	     price. "Pre-Hurricane" is an indicator variable for whether a station ///
		 lies in an area impacted by a hurricane or coastal area five ///
		 days before a hurricane warning was issued. "Hurricane" and ///
		 "Post-Hurricane" are similar indicator variables for days during a ///
		 hurricane warning and the five-days after a hurricane warning, ///
		 respectively. Standard errors are clustered at the ///
		 county for retail regressions and wholesale rack for wholesale ///
		 regressions. *, **, and ***   denote significance at ///
		 the 10\%, 5\%, and 1\% level.) ///
	coef(pre_hur "Pre-Hurricane" hur "Hurricane" post_hur "Post-Hurricane") ///
	scalars("n_stats Stations/Racks" "state_fe State FE"  ///
	        "station_fe Station/Rack FE" "year_fe Year FE" ///
			"month_fe Month-of-Year FE" "dow_fe Day-of-Week FE")  ///
			 sfmt(%8.0f) mlabels((1) (2) (3) (4)) collabels(none) 


*****************
* PANEL B AVERAGE HURRICANE IMPACTS ON RETAIL AND WHOLESALE MARGINS
*****************
eststo clear 		   
*Column 1 - Retail 1
reghdfe retail pre_hur hur post_hur wholesale CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(state_name year month dow) cluster(county_FIPS)	   
	gdistinct station_id if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "Yes"
	estadd local station_fe "No"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"
	estadd local dow_fe "Yes"	
est store A
			   
*Column 2 - Retail 2
reghdfe retail pre_hur hur post_hur wholesale CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(station_id year month dow) cluster(county_FIPS)			   
	gdistinct station_id if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "No"
	estadd local station_fe "Yes"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"	
	estadd local dow_fe "Yes"	
est store B

*Column 3 - Wholesale 1
reghdfe wholesale pre_hur hur post_hur bulk CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(state_name year month) cluster(nearRack1)
	gdistinct nearRack1 if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "Yes"
	estadd local station_fe "No"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"
	estadd local dow_fe "Yes"	
est store C
			   
*Column 4 - Wholesale 2
reghdfe wholesale pre_hur hur post_hur bulk CTST CHRT CTSN CHRN temp temp2, ///
			   absorb(nearRack1 year month) cluster(nearRack1) 		
	gdistinct nearRack1 if e(sample)
	estadd scalar n_stats = r(ndistinct) 
	estadd local state_fe "No"
	estadd local station_fe "Yes"
	estadd local year_fe "Yes"
	estadd local month_fe "Yes"		
	estadd local dow_fe "Yes"	
est store D

esttab A B C D using "$output/twfe_b_$outputdate.csv", replace label ///
     b(a2) nonumber ///
    starlevels(* 0.10 ** 0.05 *** 0.01) ///
	title(Average Effect of Hurricanes on Retail Prices, Wholesale Prices, and Margins  \label{tab:price_main1b})  ///
	cells(b(fmt(3) star) se(fmt(3) par)) ///
	drop(_cons CTST CHRT CTSN CHRN temp temp2) ///
	note(Notes: The dependent variable is station-level retail/wholesale ///
	     price. "Pre-Hurricane" is an indicator variable for whether a station ///
		 lies in an area impacted by a hurricane or coastal area five ///
		 days before a hurricane warning was issued. "Hurricane" and ///
		 "Post-Hurricane" are similar indicator variables for days during a ///
		 hurricane warning and the five-days after a hurricane warning, ///
		 respectively. Standard errors are clustered at the ///
		 county for retail regressions and wholesale rack for wholesale ///
		 regressions. *, **, and ***   denote significance at ///
		 the 10\%, 5\%, and 1\% level.) ///
	coef(pre_hur "Pre-Hurricane" hur "Hurricane" post_hur "Post-Hurricane") ///
	scalars("n_stats Stations/Racks" "state_fe State FE"  ///
	        "station_fe Station/Rack FE" "year_fe Year FE" ///
			"month_fe Month-of-Year FE" "dow_fe Day-of-Week FE")  ///
			 sfmt(%8.0f) mlabels((1) (2) (3) (4)) collabels(none) 
			 
frame change default
	  